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I. INTRODUCTION 



A. BACKGROUND 

Since September 1992, the Department of the Navy’s strategic planning document 
[Ref. 1] has officially required U.S. forces to be ready to operate in “littoral” or shallow 
water environments instead of the open-ocean climates it had based it’s warfighting 
strategies on previously. This change in policy significantly increased the potential for 
relatively inexpensive weapons, such as underwater mines, to be purchased and layed in 
the “littoral” environments U.S. Navy warships now routinely patrol. The underwater 
explosion from one of these weapons can produce shock waves or pressure pulses that, 
when applied to the large area of the ship’s hull, can have a serious negative effect on the 
ship’s structure and equipment, as well as causing serious personnel casualties. 
However, this is not a new concern. As a defensive measure, the U.S. Navy, since the 
Second World War, has been developing guidelines and specifications for the shock 
testing and hardening of shipboard equipment and systems. NAVSEA 0908-LP-000- 
3010A [Ref. 2] and MIL-S-901D [Ref. 3] are examples of this guidance. The shock 
resistance validation is then conducted through shock trials as required in OPNAVINST 
9072.2 [Ref. 4]. At this time, shock trials are the only means of testing a ship and it’s 
mission-critical systems under combat-like conditions short of an actual conflict. These 
trials are required for the lead ship of each new construction shock hardened ship class. 
These requirements have since proved their worth, when in 1991 during Operation Desert 
Storm, the USS Princeton (CG-59) struck a floating mine near the bow causing severe 
hull girder damage near the stem of the ship and only minor personnel casualties. There 
may have been massive loss of life as well as the complete loss of the ship had it not been 
for the heroic damage control efforts, as well as the “shock hardening” design 
requirements mandated by the Department of the Navy. 

Unfortunately, conducting ship shock trials can be time consuming and expensive. 
Initial planning stages for the shock trial of the USS John Paul Jones (DDG 53) began 
four years prior to the actual test date while the ship was still undergoing construction at 
Bath Iron Works, Maine. The entire endeavor involved over 50 government agencies. 
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Originally scheduled for February 1994, the shock trial was delayed 3 months due to a 
lawsuit filed against the Navy by environmentalist groups concerned over the well being 
of sea life in the testing area. When testing occurred in June 1994, only two of the four 
required tests could be carried out because of inclement weather. The remaining two 
shock tests could not be performed. Further modifications to the ship’s schedule to 
accommodate the two tests were not feasible since the three-month delay had already 
affected the ship’s post trial delivery date and deployment preparations [Ref. 5]. 

With the advent and ongoing advances in computer technology, finite element 
modeling and simulation has become a viable, less costly alternative to live fire testing. 
Finite element modeling using codes such as TrueGr/d® [Ref. 6] and MSC/PATRAN 
[Ref. 7] have enabled the generation of detailed finite element models in a timely manner. 
To analyze finite element models, highly detailed model meshes are required to provide 
the most reliable results. An important aspect of this model detail is the inclusion of the 
surrounding fluid in order to accurately capture the response of the ship caused by the 
impact of the shock wave. The fluid mesh must be constructed to mate exactly with the 
finite element mesh of the structure and must be of sufficient size to capture the bulk 
cavitation zone. However, depending on the size of the charge producing the underwater 
explosion event, as well the depth of the charge, this bulk cavitation zone can become 
quite large. The larger the cavitation zone, the greater the number of finite elements 
required to accurately model the fluid. Subsequently, greater computational memory and 
time are required to perform the structural analysis. 



B. SCOPE OF RESEARCH 

This paper investigates the effects of reducing the fluid mesh size on the accuracy 
of the structural response of a finite element model subjected to an underwater explosion. 
The model to be considered in this study is of the Navy’s Floating Shock Platform (FSP). 
The FSP is used in the acceptance of mission-essential items for installation aboard shock 
hardened ships if the size and weight of the item permits such testing [Ref. 4], Analysis 
of the model response is conducted using the LS-DYNA/USA (Underwater Shock 
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Analysis) coupled computer code [Ref.’s 13 and 14]. The purpose of this thesis is to 
compare the responses of truncated or reduced fluid mesh size models to that of a full 
size model, which is comprised of fluid volume elements within the bulk cavitation zone, 
as well as elements outside of the cavitation zone, closer to the location of the explosive 
charge. With reliable results, computer simulation of ship shock trials could become a 
dependable, cost effective, and time efficient manner for validating surface ship shock 
hardening requirements. 
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II. UNDERWATER EXPLOSIONS 



A. UNDERWATER SHOCK PHENOMENA 

An underwater explosion occurs in a complex sequence of events. When a high 
explosive, such as TNT or HBX-1 is detonated, the original solid material of the 
explosive is converted into a very high temperature and pressure gas within nanoseconds 
(on the order of 3000°C and 50000 atm.) [Ref. 8]. The pressure wave that is formed 
originates in one section of the explosive and propagates throughout the remainder of the 
explosive. As this pressure wave propagates, it initiates the chemical reaction that creates 
more pressure waves. The pressure wave velocity steadily increases within the solid 
explosive until it exceeds the speed of sound in the explosive, creating a shock wave. 
The shock wave propagates through the solid at a constant speed and then, with the high 
temperature and pressure behind the shock front, into the surrounding medium [Ref. 9]. 

The high-pressure gas that results from the explosion rapidly expands outward in 
a radial manner, Fig. 1, and imparts an outward velocity on the surrounding water as 




Figure 1. Gas Bubble and Shock Wave from an Underwater Explosion 

well. Initially, the pressure is much greater than the atmospheric and hydrostatic pressure 
that opposes it and is therefore compressive in nature. At detonation, the pressure rise 
produces a steep fronted discontinuous wave, which decays exponentially with time as 
shown in Fig. 2. Duration of the pressure disturbance lasts only a few milliseconds. The 
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shock wave is assumed to propagate at several times that of the speed of sound in water, 
approximately 5,000 ft/sec, near the charge, which then falls rapidly to acoustic velocity 
as it travels outward in the water. 



34,000 psi 




5' 45' 50' 495’ 500’ 



Figure 2. Shock Wave Profiles From a 300 lb. TNT Charge [Ref. 8] 



Additionally, the pressure profile of the shock wave is proportional to the inverse of the 
distance from the charge, 1/d, and the w'ave profile gradually broadens as it spreads out 
[Ref. 8]. Empirical equations have been determined to define the profile of the shock 
wave. These relations enable calculation of the pressure profile of the shock wave (P(t)), 
the maximum pressure of the wave (P m ax), the shock wave decay constant (0), the bubble 
period (T), and the maximum bubble radius (A max ). 

P(t) = P„„e‘^ (Ps>) (2-1) 



P max = K, 

max 1 
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A„ = K« — | (ft) ( 2 . 5 ) 

(D + 33)1 

Other variables in the equations are: 

W = Charge weight (lbf) 

R = Standoff distance (ft) 

D = Charge depth (ft) 

ti = arrival time of shock wave (msec) 

t = time of interest (msec) 

Ki, K2, K5, Kg, Ai, A2 = Shock wave parameters 

Through calculation, it can be determined that P max decreases by approximately one-third 
after one decay constant. 

Subsequent pressure waves or bubble pulses are generated by the oscillation of 
the gas bubble created by the underwater explosion. The peak pressure of the first bubble 
pulse is approximately 10-20% of the shock wave, but is of greater duration making the 
area under both pressure curves similar [Ref. 9 ]. The bubble expands until dynamic 
equilibrium is reached . The bubble then contracts until dynamic equilibrium is again 
reached, followed by another expansion. This oscillation sequence continues until the 
energy of the reaction is dissipated or the bubble reaches the free surface or impacts the 
target. 

Depending on the charge location relative to the surface and the bottom, other 
effects are characteristic of an underwater shock. Bottom reflection is the bouncing of 
the shock wave off of the bottom of the body of water; a compressive wave. Refraction 
causes the shock wave to travel through the bottom of the body of water before emerging 
again; also a compressive wave. In reasonably deep water, these two effects are not 
usually an issue for surface vessels. 
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Free surface reflection is a very important effect, however. This reflected wave is 
tensile in nature, as opposed to the other compressive wave effects, and is produced from 
the rarefaction of the shock wave from the free surface. This rarefaction wave 
contributes to the creation of bulk cavitation. 



B. FLUID-STRUCTURE INTERACTION 

When an object such as a ship or submarine is in the vicinity of an underwater 
explosion, the shock pressure pulses produced by the explosion impinge upon the surface 
of the structure. A fluid-structure interaction takes place as the pressure pulse acts upon 
the flexible surface of the structure. This dynamic response of the a linear elastic 
structure in the fluid can be expressed by: 

[M s ]{x} + [CJ{x} + [KJ{x} = {f} (2.6) 

where [M s ] is the mass matrix, [C s ] is the dampening matrix, [K s ] is the stiffness matrix, 
{x} is the acceleration vector, {x} is the velocity vector, and {x} is the displacement 
vector of the structure and {f} is the external force vector. In the case of a submerged 
structure excited by an acoustic wave, {f} is given by: 

{f} = -[G][A f ]({ Pl } + {p s }) + {f D } (2.7) 

where [G] is the transformation matrix that relates the surface nodal forces of the fluid 
and structure, [Af] is the diagonal area matrix associated with the fluid elements, 
{pi}=incident wave nodal pressure vector, and {ps} -scattered wave nodal pressure vector 
[Ref. 11]. 

The Doubly Asymptotic Approximation (DAA) method was developed in 1971 
[Ref. 10], which described a matrix of differential equations in time for the 
approximation of acoustic fluid-structure interaction. This approximation was so named 
because it is accurate at both low and high frequencies and at early and late times. The 
DAA represents the surrounding fluid of the structure through the interaction of state 
variables pertaining only to the structure’s wet surface. The fluid equation of the DAA is 
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[M f ](P S } + pc[A f ]{p s } = pc[M f ]{u s } (2.8) 

where [Mf] is the symmetric fluid mass matrix for the wet-surface fluid mesh, {ps} and 
{p s } are the nodal pressure vector and its first time derivative of the scattered wave, c is 
the acoustic velocity of water, [Af] is the diagonal area matrix associated with the fluid 
elements, and {u s } is the scattered wave velocity vector. An object in the fluid will have 
a structural response defined as follows in Eq. (2.9), 

[M s ]{x} + [CJ{x} + [K s ]{x} = {f} (2.9) 

where [M s ] is the mass matrix, [C s ] is the dampening matrix, [K s ] is the stiffness matrix, 
{x} is the acceleration vector, {x} is the velocity vector, and {x} is the displacement 
vector of the structure and {f} is the external force vector. 

The kinematic compatibility relation can then be applied to relate (us) to the 
structural response, 

[G] T {x} = (u,} + {u s } (2.10) 

The “T” superscript indicates the transpose of the matrix. This equation is an expression 
of the constraint that the normal fluid particle velocity must match the normal structural 
velocity on the structure wetted surface. 

Substituting Equation (2.7) into (2.6) and Equation (2.10) can be substituted into 
(2.8) resulting in Equations (2.1 1) and (2.12), 

[MJ{x} + [CJ{x} + [KJ{x} = -[G][A,]({ Pl } + {p s }) (2.11) 

and 

[M f ] {p s } + pc[A f ] ( Ps } = pc[M f ]([G] T {x} - {u , } ) (2.12) 

The Underwater Shock Analysis (USA) code solves Eqs. (2.11) and (2.12) 
simultaneously by using a staggered solution procedure that is unconditionally stable 
with respect to the time step used [Ref. 11]. Once this system of equations is solved, 
desired response results such as displacement, velocity, and acceleration can be studied. 
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C. CAVITATION 

Two types of cavitation can occur during an UNDEX event. "Local cavitation" 
occurs at the fluid-structure interface and "bulk cavitation" occurs near the free surface 
and can cover a relatively large area. A discussion of both forms of cavitation follow 
below. 

1. Local Cavitation 

Taylor flat plate theory, the simplest case of fluid-structure interaction is used to 
illustrate how local cavitation occurs. In this case, an infinite, air-backed plate is acted 
upon by an incident plane shock wave as shown below in Fig. 3. 




Shock Wave Approaching Shock Wave Reflecting 



Figure 3. Taylor Plate Subjected to a Plane Wave [Ref. 9] 



Once the shock wave strikes the plate, a reflected shock wave leaves the plate. 
According to Newton’s second law of motion as shown in Eq. (2.13), 

m^ = P,+P, (2.13) 

at 

where m is the mass of the plate per unit area, u is the velocity of the plate after being 
subjected to the shock wave, Pi(t) is the incident wave pressure and P 2 (t) is the reflected, 
or scattered, wave pressure. Define the fluid particle velocities behind the incident and 
reflected shock waves as Ui(t) and U2(t). The velocity of the plate is then defined by Eq. 
(2.14), 
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u(t) = u,(t)-u 2 (t) (2.14) 

For the one dimensional plane wave, the wave equation is P = pCu. It follows that the 
incident and reflected shock wave pressures, Eqs. (2.15) and (2.16), are 

P. =pCu, (2.15) 

P 2 =pCu 2 (2.16) 

where p is the fluid density and C is the acoustic velocity in water. Substituting the 
above pressure Eqs (2.15) and (2.16) into the velocity Eq. (2.14) results in the incident 
shock pressure being defined as Eq. (2.17) 

P i( t ) = P«,».e'“ (2.17) 

Solving for the reflected shock pressure yields Eq. (2.18), 

P 2 (t) = P,-pCu = P max e 9 -pCu (2.18) 

where t is the time after the shock wave arrives at the target. Now the equation of 
motion, Eq. (2.13) can be rewritten as Eq. (2.19), 



m (^) + PCU = (2 ' 19) 

which is a first order, linear differential equation. The solution, u(t), of the differential 
equation is expressed in Eq. (2.20) as 



2P e -i 

u = _ £f™* u .. r e e - e e l 

m(l-p) L 



( 2 . 20 ) 



with P = pC9/m and t>0. The total pressure that impinges on the plate is defined as Eq. 

( 2 . 21 ), 



P +P =P I 

1 1 ~ 1 2 1 max I 



1-P 



2(3 

l-(3 



.Pi 

e e ] 



( 2 . 21 ) 
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As the value of p becomes larger, as in the case of a lightweight plate, the total pressure 
will become negative at a very early time. However, since water cannot support tension, 
negative pressure cannot exist. Therefore, as the water pressure reduces to vapor 
pressure at the surface of the plate, cavitation occurs. At this point, the pressure in front 
of the plate has been cut off and the plate has reached its maximum velocity [Ref. 9] . 

A ship's hull can be easily generalized as a Taylor flat plate. Local cavitation is 
likely to occur along the hull where the pressure pulse from the UNDEX impinges with 
sufficient force and the hull plating P value is large enough to make the net pressure 
negative. 

2. Bulk Cavitation 

The incident shock wave is compressive in nature. A rarefaction wave, which is 
tensile in nature, is created when the shock wave is reflected from the free surface. Since 
water cannot sustain a significant amount of tension, cavitation will occur when the 
pressure drops to zero or below. Upon cavitation, the water pressure rises to the vapor 
pressure of water, approximately 0.3 psi. This cavitated region created by the rarefaction 
wave is known as the bulk cavitation zone. It consists of an upper and lower boundary 
and its extent is dependent on the charge size, type, and depth [Ref.’s 18 - 20]. 

Figure 4 shows a typical bulk cavitation zone. The cavitation zone is symmetric 
about the y-axis in the figure; typically only one-half is shown due to the symmetry. The 
water particles behind the shock wave front at the time of cavitation have velocities 
depending on their location relative to the charge and the free surface. Water particles 
near the free surface, for example, will have a primarily vertical velocity at cavitation. 
As the reflected wave passes, the particles will be acted upon by gravity and atmospheric 
pressure. 

The upper cavitation boundary is the set of points where the rarefaction wave 
passes and reduces the absolute pressure to zero or a negative value. The region will 
remain cavitated as long as the pressure remains below the vapor pressure. The total or 
absolute pressure, which determines the upper boundary, is a combination of atmospheric 
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pressure, hydrostatic pressure, incident shock wave pressure, and rarefaction wave 
pressure. 




The lower cavitation boundary is determined by equating the decay rate of the 
breaking pressure to the decay rate of the total absolute pressure. The breaking pressure 
is the rarefaction wave pressure that reduces a particular location of a fluid to the point of 
cavitation pressure, or zero psi. 

The upper and lower cavitation boundaries can be calculated from Equations 
(2.22) and (2.23), respectively [Ref. 19]. Any point which satisfies F(x,y) and G(x,y) = 0 
determines the bulk cavitation boundary. 
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The variables in Equations (2.22) and (2.23) are: 



x, y = horizontal range and vertical depth of the point 

ri = standoff distance from the charge to the point 

12 = standoff distance from the image charge to the point 

C = acoustic velocity in the water 

D = charge depth 

9 = decay constant 

y = weight density of water 

Pa = atmospheric pressure 

W = charge weight 

P, = P(t), Equation (2. 1 ) 

9 = Equation (2.3) 

Kj, A] = shock wave parameters 

Figure 5 shows the charge geometry for the above two equations. 




Figure 5. Charge Geometry for Bulk Cavitation Equations [Ref. 9] 

Appendix A provides a MATLAB m-file [Ref. 12] that calculates and plots the 
bulk cavitation zone for a user supplied charge weight (of HBX-1) and depth by solving 
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Equations (2.22) and (2.23). Figure 6 provides an example of cavitation curves generated 
using the program for two different charge weights at three different depths. 
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Figure 6. Bulk Cavitation Zones for HBX-1 Charges at the Following Depths: 

- 100 ft, - 200 ft, 300 ft 
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III. MODELING 



Modeling and testing involves model construction and pre-processing, analysis 
and solution, and post-processing programs. A flow chart of the model building and 
testing procedure is shown below in Figure 7. 



Model Construction 
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Figure 7. Flow Chart. Model Construction and Testing 
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A. MODEL CONSTRUCTION AND PRE-PROCESSING 



1. Floating Shock Platform (FSP) Structural Model 

The model of the Navy Standard Floating Shock Platform used for the 
heavyweight shock simulation was constructed using a finite element mesh generation 
program called Tru zGrid [Ref. 6]. 

The model was developed using information detailed in the Military Specification 
for Shock Tests [Ref. 3]. The dimensions of the Navy Standard Floating Shock Platform 
are 288-in long, 192 -in wide, 72-in deep on outside, and 36-in deep on the interior. For 
this application, it was possible to model only half of the platform due to symmetry. This 
was done to significantly reduce the computation time. This made the final model 
dimensions 144-in long, 192-in wide, 72-in deep on outside, and 36-in deep on the 
interior. The shell plating for the exterior of the modeled using thickness of l A and 1-in 
steel having weight densities of 0.283, 0.317, and 0.340 lbf/in 3 , Young’s Modulus of 
30x1 0 6 psi, and a Poisson’s ratio of 0.3. The overall finite element mesh consists of 3669 
nodes and 4116 quadrilateral (shell) elements. Table 1, Figures 8 and 9 show the model 
particulars. Figure 10 shows the overall finite element model. 



Length 


144-in 


Beam 


192-in 


Depth (inside) 


36-in 


Depth (outside) 


72-in 


Design Waterline 


36-in 


Plating/Stiffener Material 


40# S.T.S. and 20.4# HY-80 Steel 


Plating/Stiffener Thickness 


!4-in and 1-in 


No. of Nodes 


3669 


No. of Shell Elements 


4116 



Table 1. Model Specifications 
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Tod View 



Figure 8. Model Specifications 
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Side View 



Figure 9. Model Specifications 




Figure 10. Finite Element Mesh 
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2. Fluid Modeling 

Upon completion of the structural model, the next step was to design the fluid 
mesh. Element extrusion was performed using TmeGricTs [Ref. 6] extrude command. 
Appendix D describes the use of this function in detail. The fluid mesh consists of 8- 
noded solid elements. LS-DYNA’s Material Type 90 (acoustic pressure element) is used 
to model the pressure wave transmission properties of water [Ref 13]. Three different 
fluid size models were constructed for analysis. The extent (in the x and z directions) of 
the first fluid model was 672-in by 372-in, respectively. The depth of the mesh (under 
the FSP) was set to 108-in to ensure that the lower fluid boundary extended beyond the 
computed bulk cavitation zone (to be discussed later). This mesh (referred to as the full 
fluid size model) contains 57768 8-noded elements and 65808 total nodes and is shown in 
Figure 11. The size and complexity of a fluid mesh of this size requires extensive 
computational power to run a shock simulation. It is therefore beneficial to be able to 
reduce the size of this fluid mesh to save computational time. 

The second fluid size model is a truncated or reduced size fluid mesh. The extent 
(in the x and z directions) of this mesh was set to 552-in by 276-in, respectively. The 
depth of this mesh (under the FSP) was set to 72-in, which places approximately half of 
the fluid mesh lower boundary outside the bulk cavitation zone. This mesh contains 
47256 8-noded elements and 54636 total nodes and is shown in Figure 12. This fluid 
model will be referred to as the first fluid truncation. 

The third fluid size model is a further truncation of the first fluid size model. It’s 
extent (in the x and z directions) was set to 432-in by 228-in, respectively. The depth of 
this mesh (under the FSP) was set to 36-in, which places the entire fluid mesh within the 
bulk cavitation zone. This mesh contains 21010 8-noded elements and 27084 total nodes 
and is shown in Figure 13. This fluid model will be referred to as the second fluid 
truncation. A summary of each fluid size model is shown in Table 2. 

An important aspect of a fluid mesh is the element size next to the structural 
mesh. For the cavitation analysis using the USA code, the critical element size is 
determined by the following equation [Ref. 15]: 
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2Ep< 5 

PA 



Where p = density of water, D = thickness of the fluid element in the direction normal to 
the wetted surface of the structure, p s = density of the submerged structure, and t s = 
thickness of the submerged structure. It can be shown for this model that the critical 
element thickness, D, is 10 inches. The elements adjacent to this structural model were 
set a value less than this value of thickness. 





Full Fluid Size Model 


First Truncation 


Second Truncation 


Length (x-dir) 


372-in 


276-in 


228-in 


Width (z-dir) 


672-in 


552-in 


432-in 


Depth (y-dir) 


108-in 


72-in 


36-in 


Number of Nodes 


65808 


54636 


27084 


Number of Elements 


57768 


47256 


21010 



Table 2. Summary of Finite Element Fluid Models 
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Figure 11. Full Fluid Size Model 
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Figure 12. First Fluid Truncation 
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Figure 13. Second Fluid Truncation 
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B. ANALYSIS AND SOLUTION 



1. Analysis Program Description 

The finite element model must be translated into LS-DYNA keyword format in 
order to perform the analysis, since the LS-DYNA/USA coupled code is used. The USA 
code performs the majority of the work (formulation of the fluid-structure interaction 
matrices) and LS-DYNA is used to perform the time integration solution for the structure. 
LS-DYNA is a non-linear three-dimensional structural analysis code [Ref. 13]. The USA 
code consists of three main modules: FLUMAS, AUGMAT, and TIMINT. 

FLUMAS is the first USA module required to be run. FLUMAS generates the 
fluid mass matrix for the submerged portion of the structure [Ref. 14]. The fluid mesh 
data, as well as the transformation coefficients that relate both the structural and fluid 
degrees of freedom on the wetted surface are generated, including the nodal weights for 
the fluid element pressure forces and the direction cosines for the normal pressure force. 
The fluid area matrix is diagonal and the fluid mass matrix is fully symmetric. 

AUGMAT is the second module to be run. This module takes the data generated 
by the FLUMAS and initial LS-DYNA runs to construct specific constants and arrays 
utilized in the staggered solution procedure for the actual transient response analysis 
[Ref. 14]. The augmented interaction equations are formed from Equations (2.11) and 
(2.12). These two equations may be solved simultaneously at each time step, but this 
solution method can be very computationally expensive. The USA code uses a staggered 
solution procedure to achieve an efficient solution. 

TIMINT performs the direct numerical time integration and also handles the 
computation of the UNDEX parameters, such as the shock wave pressure profile. The 
structural and fluid response equations are solved separately at each time step through the 
extrapolation of the coupling terms for the two systems. LS-DYNA is used to solve the 
structural equations and the TIMINT run solves the fluid equations. A result of using the 
staggered solution procedure mentioned previously, is that LS-DYNA and TIMINT can 
each have a different time step assigned. Although, the general practice is to use the 
same time step in both computations. Despite using an unconditionally stable solution 
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scheme, the TIMINT time step must be set small enough to accurately capture the fluid 
system response. Additionally, it should be mentioned that LS-DYNA uses a central- 
difference integration method, which is conditionally stable. The LS-DYNA time step 
must be less than or equal to the critical time step for the structural finite element mesh or 
numerical instability will result. Overall, this step of the solution procedure is the most 
time consuming and computationally expensive. 

Appendix B provides example input decks for each of the three USA modules as 
well as an example LS-DYNA Keyword input deck. 

2. Test Description 

Two different charge location geometries were used in the shock simulation runs 
for this study. They are similar to charge locations specified for heavyweight shock 
testing of standard FSP’s by Military Specifications for Shock Tests [Ref. 3]. 

A charge consisting of 60 lb. HBX-1 was chosen because of its specification in 
the same reference. One attack geometry placed the charge offset from the side of the 
FSP at the plane of symmetry (y-z plane) by 240-in. The charge depth is 288-in, with a 
standoff distance of 375-in. Figure 14 shows this attack geometry and Table 3 shows a 
summary of the UNDEX parameters of the explosion. The bulk cavitation zone was 
computed using the MATLAB program in Appendix A and is shown in Figure 15 along 
with the three fluid size mesh model boundaries. The second attack geometry consisted 
of the same 60 lb. HBX-1 charge, but placed directly under the FSP at the plane of 
symmetry (y-z plane) at a depth of 288-in. This resulted in a standoff distance of 252-in. 
Figure 16 shows this geometry and Table 4 summarizes the UNDEX parameters for this 
configuration. Once again, the bulk cavitation zone is shown in it’s entirety and also on a 
larger scale which includes the three fluid size mesh model boundaries with respect to the 
FSP in Figure 17. 
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Y 




HBX-1 



Figure 14. Offset Charge Test Geometry 



Standoff Distance 


375.0 in 


P max 


2077 psi 


0 


0.3662 msec 


T 


0.64 sec 


A m ax 


172.60 in 



Table 3. UNDEX Parameters for Offset Charge 
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Entire Bulk Cavitation Zone 




FSP and Fluid Mesh Boundary Location in Bulk Cavitation Zone 




Figure 15. FSP in Bulk Cavitation Zone for a 60-lb. Charge (Offset) at 288-in 
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60 lb 
HBX-1 



Figure 16. Charge Under FSP Test Geometry 



Standoff Distance 


252.0 in 


Pmax 


3271 psi 


e 


0.3320 msec 


T 


0.64 sec 


A m ax 


172.60 in 



Table 4. UNDEX Parameters for Charge Under FSP 
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Feet . Feet 



Entire Bulk Cavitation Zone 





Figure 17. FSP in Bulk Cavitation Zone for a 20-lb Charge (Under) at 288-in 
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C. POSTPROCESSING 

The solution data is output into two formats for analysis: binary and ASCII. The 
binary data files created by the LS-DYNA/USA runs contain the model’s finite element 
response information. LS-TAURUS [Ref 13] and Glview [Ref. 16] can both be used for 
three-dimensional response visualization. They provide powerful animation and image 
generation features, including the display of velocity, acceleration, and element pressure 
data. Additionally, both post-processing programs enable the user to observe the shock 
wave propagation through the fluid to the structure. LS-TAURUS has the capability of 
extracting ASCII solution data and writing it to a separate ASCII file for later evaluation. 
Appendix E provides some useful LS-TAURUS commands for model post-processing. 

The ASCII data extracted by the LS-TAURUS post-processing program can be 
plotted and manipulated using UERD (Underwater Explosion Research Division) Tool 
software. This program is a PC-based plotting tool, which plots ASCII input files, 
provides standard graphing functions, as well as providing a variety of data manipulation 
features. These features include, but are not limited to curve integration and derivation of 
shock spectra. 
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IV. FSP SHOCK SIMULATION RESULTS 



Three node points were used for comparison between the different fluid sized 
models in the simulation. The velocity response was analyzed at these nodes, specifically 
in the vertical (y) direction. These nodes are all located on the interior (floor) of the FSP. 
One node at the middle of the platform (closest to the symmetry plane), one at the interior 
end adjacent to the aft bulkhead, and one node V 2 the distance between the other two 
nodes were chosen to provide a response at positions varying in stiffness. These nodes 
and their ID numbers are shown in Figures 8 and 10. 

A. CHARGE OFFSET FROM FSP 

The first case examined was the geometry with the charge offset from the FSP. It 
was possible to run the simulation by placing the DAA boundary directly on the wetted 
surface of the structure, DAA boundary on the fluid mesh with cavitation flag off, or 
DAA boundary on the fluid mesh with the cavitation flag on. Results from a previous 
thesis [Ref. 17] concluded that the method of approximation the response by placing the 
DAA boundary on the fluid mesh with cavitation flag on is the most accurate. Therefore, 
this was the method chosen for this simulation. The simulation was run for all three fluid 
size mesh models in order to compare the effects of a truncated mesh on the response of 
the FSP. 

The responses are plotted using the full fluid mesh size model results as a 
reference. Truncated fluid mesh size model responses are then plotted one at a time 
along with the full fluid mesh size response at each node for easier and clearer 
comparison. All simulations were run on an SGI Octane with two 195 MHz processors, 
1 .344 Gigabytes of RAM, and 23 Gigabytes of hard drive storage capacity. LS-DYNA 
version 940.2b and USA+ version 980 IE were the simulation codes used. 

For the full fluid mesh size model, a total of 3 hours was required for FLUMAS 
module to complete its computations. The AUGMAT module took 1 hour, 45 minutes to 
run, and the TIMINT module took 3 hours, 55 minutes to complete its calculations (with 
a time step of 2 x 10‘ 5 seconds for 80 msec of data) for a total run time of 8 hours, 50 



minutes. There was a considerable savings in computation time experienced for each 
truncated fluid mesh model. For the second fluid mesh truncated model, a total of 39 
minutes was required for the FLUMAS module to complete its computations. The 
AUGMAT module took only 26 minutes to run, and the TIMINT module took just 1 
hour, 29 minutes to complete its calculations (with the same time step and data duration) 
for a total run time of 2 hours, 34 minutes. The computation time was V 3 of the time 
necessary to run the full fluid mesh size model. 

The vertical velocity response of the structure is what was expected. The velocity 
increases rapidly (approximately 5 to 6 msec) to a peak value and then rapidly decreases 
and increases in a series of steps until it settles down. The response does not settle out at 
a value of zero but instead slightly above that value due to rigid body motion of the 
structure. The reason for this behavior is from the fact that the incident shock wave 
impacts the structure with a very high pressure (close to 2100 psi) at approximately 5 to 6 
msec from charge detonation and forces the structure rapidly upward. The structure is 
then quickly pulled downward as the shock wave reaches the free surface and a tensile 
wave is generated. This wave causes the pressure to decrease rapidly to zero psi, and 
cavitation occurs. Once cavitation occurs, the FSP is allowed to “break free” of the fluid, 
and the velocity of the structure again increases, albeit to a lower magnitude until the 
cavitating fluid contacts the structure pulling it downward again. This cycle continues 
past the point where the velocity of the structure is even negative momentarily to a value 
of just above zero, once again due to rigid body motion of the FSP. 

The effect of fluid mesh truncation on the response of the structure is very slight. 
This can be seen in Figures 18 and 19. There is only a slight divergence after 
approximately 40 msec. Even with a fluid truncated model of approximately Vg the size, 
the response is extremely accurate. 

The athwartships (z-direction) velocity was also examined at the three structural 
nodes. Once again, there is fairly good correlation between the full fluid size meshed 
model and the first and second fluid truncated models. 
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The pressure at the middle and bottom of the fluid mesh underneath the FSP was 
examined. These pressure plots are included in Figures 22 and 23 and are element 
pressures taken directly below the structure. Cavitation can be seen to occur almost 
immediately underneath the model. Cavitation occurs almost continuously until 
approximately 40 msec, when the effects of the initial shock wave and rarefaction waves 
have subsided allowing the fluid pressure to again reach a state above zero. This 
response correlates well with the velocity response seen in the FSP in Figures 18 through 
21. Additionally, pressure profiles at the depths examined were expected to show 
cavitation for a large duration (after the initial shock wave passed) since these depths fall 
within the bulk cavitation zone computed and shown in Figure 15. 

It should be pointed out that the magnitude of the initial shock wave pressure 
shown in the pressure profile plots in Figures 22 and 23 was far below what was expected 
(approximately half of the expected value). This can be explained by the fact that the 
decay constant of the incident shock wave has a value of 0.3662 msec as shown in Table 
3. With a decay constant this high, the pressure wave will decay over a distance of 
approximately 21.75 inches (assuming a fluid sound speed of 4952.4 ft/sec). The fluid 
mesh refinement used for this study has a diagonal distance across elements of 14 inches. 
This distance allows the incident shock wave pressure to decay to value less than half of 
its actual pressure before it is imparted on the element adjacent to it. This is expected to 
cause an deviance in fluid pressure of approximately 50 percent when the value is 
extracted from the LS-DYNA/USA program run data files. Indeed, this is what is seen in 
Figures 22 and 23. 

This difference between program calculated values of pressure and actual 
expected values was not considered a major concern, since the majority of the FSP’s 
response is developed from the fluid particle’s impulse imparted on the structure and not 
the peak pressure of the incident wave on it. Also, for this study, the interest was more 
toward the comparison of the FSP’s structural response and fluid pressure profile when a 
fluid mesh of a given size was truncated. 
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A solution to resolve the fluid pressure decay across the fluid mesh element is to 
reduce the size of the element used. For this case, when the decay constant is high, the 
fluid element dimensions must be reduced to a size, which does not allow a significant 
reduction in peak pressure across the element during program integration. 
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Figure 18. Vertical Velocity Response w/Charge Offset from FSP 
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Figure 19. Vertical Velocity Response w/Charge Offset from FSP 
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Figure 20. Athwartships (Z-Dir) Velocity Response w/Charge Offset from FSP 



39 



Node 1576 

Athwart ships Velocity (Z-Dir) 




Time (msec) 

■ Full Fluid Size Model — — Second Fluid Truncation 



Node 2445 

A th wardships Velocity <Z-E>ir) 




Time (msec) 

' Full Fluid Size Model Second Fluid Truncation 

Node 2471 



A thwart ships Velocity (Z-Dir) 




Time (msec) 

— Full Fluid Size Model Second Fluid Truncation 



Figure 21. Athwartships (Z-Dir) Velocity Response w/Charge Offset from FSP 
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Figure 22. Fluid Mesh Pressure Profiles w/Charge Offset from FSP 
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Figure 23. Fluid Mesh Pressure Profiles w/Charge Offset from FSP 
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B. CHARGE UNDER FSP 



The shock simulations were also run and compared for the geometry with the 
charge placed directly under the FSP along the plane of symmetry (y-z plane). The same 
combination of runs was conducted as for the charge offset from the FSP. The resulting 
responses are quite similar in nature to those described in section A of this chapter. 

Vertical velocity response profiles are plotted in Figures 24 and 25. Once again, 
the velocity response increases rapidly to a peak value and then rapidly decreases and 
increases in a series of steps until it settles down. The response for this case takes 
slightly longer to reach a somewhat steady state because the magnitude of the incident 
pressure wave is a larger value (3271, vice 2077 psi). However, as in the offset charge 
geometry, the structural velocity settles out to a value above zero because of the rigid 
body motion of the FSP. 

The effect of the fluid mesh truncation on the response of the structure is slight as 
in the previous geometry. The slight divergence encountered occurs after approximately 
60 msec. However, even with the extent of truncation between the full fluid size model 
and the second fluid truncated model, the deviance is not large. 

Again, pressure at the middle and bottom the fluid mesh underneath the FSP was 
compared. These pressure profile plots are shown in Figures 26 and 27. Cavitation 
occurs as expected, almost immediately (following the initial shock wave pressure pulse). 
In this geometry, it can be more readily seen that the pressure deeper in the fluid mesh 
experiences cavitation for a much shorter duration of time. In fact, cavitation subsides 
within approximately 20 msec at depth of 29-in. as seen in Figure 27. Cavitation does 
indeed occur, as predicted, since this depth falls within the computed bulk cavitation zone 
plotted in Figure 17. The bulk cavitation zone depth directly under the FSP was 
computed to a depth of approximately 72 inches. Therefore, cavitation is expected to 
occur in each of the pressure profile comparison plots. The deeper pressure profile plot 
comparison at 62 inches below the FSP shown in Figure 26 indicates cavitation occurs 
for a much shorter duration (approximately 10 msec) as expected. 
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As was previously discussed in section A of this chapter, the magnitude of the 
initial shock wave pressure shown in Figures 26 and 27 was much less than expected. 
The cause for the deviance and solution is the same as was discussed previously. Once 
again, because the basis of this study was structural response comparison when a given 
fluid size mesh is truncated, this was not seen as a major concern, but one, which must 
still be addressed, as was done here. 
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Figure 24. Vertical Velocity Response w/Charge Under FSP 
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Figure 25. Vertical Velocity Response w/Charge Under FSP 
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Figure 26. Fluid Mesh Pressure Profiles w/Charge Under FSP 
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Figure 27. Fluid Mesh Pressure Profiles w/Charge Under FSP 
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V. CONCLUSIONS AND RECOMMENDATIONS 



This thesis investigated the effect of reducing or truncating the fluid mesh in a 
fluid-structure finite element model on the accuracy of a FSP’s response to an underwater 
explosion. There is no question that when practical, a fluid mesh size large enough to 
capture the extent of a computed bulk cavitation zone should be used to obtain the most 
accurate structural response possible. 

However, when computer extent and computational time is an issue, as is 
increasingly more probable as finite element models become more complex, it is of high 
importance to develop a means to obtain accurate results while minimizing computational 
effort. Fluid mesh truncation was explored as a means to reduce the computation time of 
LS-DYNA/USA coupled code to calculate a solution for the structural response to an 
underwater explosion event. 

The structural velocity response results of the FSP with truncated fluid size 
models agreed closely with the response obtained by modeling a larger size fluid mesh. 
The responses were not only in close agreement with the other model runs, but the values 
obtained agreed with predicted behavior based on the physics of the situation. Slight 
deviation in model results may be attributed to minor differences in program input 
parameters. In any case, fluid mesh truncation clearly demonstrated that the effects of 
fluid mesh truncation are very slight and that it is a viable option for reducing the 
computation time of solving otherwise large finite element fluid-structure models. 

It is recommended that additional studies be conducted to more fully examine the 
effects of fluid mesh truncation in underwater explosion simulations. Specifically, the 
following areas should be studied: 

1. Vary the fluid size adjacent to the FSP to determine if there is a minimum 
amount of fluid required for modeling to obtain accurate structural response. 

2. Vary the fluid element size to examine the effect on fluid element pressure 
accuracy with expected results. 
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3. Model selected equipment mounted on FSP in order to compare live fire 
testing results with computed results. 

4. Apply fluid mesh truncation to more complex models, such as U.S. Navy 
Warships, to investigate the accuracy of the model response to known live fire 
testing data. 
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APPENDIX A. BULK CAVITATION PROGRAM 



The following program code calculates the bulk cavitation zone by solving 
Equations (2.22) and (2.23). The code is written for MATLAB Version 5.2 [Ref. 12]. 

% Bulk Cavitation Program Using HBX-1 as Charge Type 
% Written by: James R. Smith 
% MATLAB Version 5.2 



% This program is used to compute both the upper and lower cavitation 
% boundaries for a given input of HBX-1 charge weight and charge depth. 
% This program can easily be modified to handle other explosive types. 

% However, shock wave parameters (Kl, K2, Al, A2) would need to be 
% updated. "Standard" atm pressure, seawater specific weight and speed 
% of sound are used. These values can be changed to achieve desired 
% accuracy for a particular problem. This program is set up as is to 
% handle cavitation boundaries out to 2200ft (x) and to a depth of 50ft 
% (y) . These parameters can be changed as required. 

% Program Start 

clear all; 



% Defining constants 



Pa = 14.7; 

C = 4900; 
gamma = 0.037031; 



% atmospheric pressure (psi) 
% acoustic velocity (ft/s) @ 
% weight density of seawater 



64. 4F ( 18C) 

0 68F and 1 atm (lbm/in A 3) 



Kl = 22347.6; 
K2 = 0.056; 

Al = 1.144; 

A2 = -0.247; 



% Shock Wave 
% Parameters 
% explosive type 
% specific 



% Prompting user to enter charge weight and depth 
disp ( ' ' ) 

W = input ( f Enter the charge weight of HBX-1 in lbf : '); 
D = input (' Enter the charge depth in feet: 1 ) ; 



W_text = num2str(W); 

D__text = num2str(D); 

disp ( 1 ') 

disp ([' Calculating the bulk cavitation zone for a f ,W_text,'lb HBX-1 
explosion 1 ] ) 

disp (['at a depth of f ,D__text,' feet 1 ]) 
disp ( 1 ' ) 



% Beginning calculation of upper and lower cavitation boundaries 

% theta = decay constant 
% x = horizontal distance 
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% y = vertical distance 

% rl,R = standoff distance from charge to point 
% r2 = standoff distance from image charge to point 

% Pi = incident shock wave pressure at tc 

data_u = []; % Creating a matrix to store upper boundary data 

data_l - []; % Creating a matrix to store lower boundary data 

% Calculating Upper Boundary 
for x = 0:600 

for y = 0:0.1:50 

rl = sqrt ( (D-y) A 2+x A 2) ; 
r2 = sqrt ( (D+y) A 2+x A 2) ; 

theta = K2*W A (1/3) * (W A (1/3) /rl) A A2/1000; 

F — (Kl* (W A (1/3) /rl) A Al*exp (- ( r2- rl ) / (C*theta) ) ) 

+Pa+ ( gamma *y* 12 ) - (Kl* (W A (1/3) /r2) A A1) ; 

if F <= 0 % Testing for cavitation 
data_u = [data_u; F x - (y) ] ; 
break, end 

end 

end 

% Calculating lower boundary 
for x = 0 : (length (data_u) -1 ) 
for y = 0:0.1:50 

rl = sqrt ( ( D-y) A 2+x A 2) ; 
r2 = sqrt ( (D+y) A 2+x A 2) ; 

theta = K2*W A (1/3) * (W A (1/3) /rl) A A2/1000; 

Pi = Kl*(W A (l/3)/rl) A Al*exp(-(r2-rl)/(C*theta)); 

G = - (Pi/ (C*theta) ) * (l+( ( (r2- 

(2*D* (D+y) /r2 ) ) /rl) * ( ( (A2*r2) /rl)-A2-l) ) )- ( (Al*Pi) /rl A 2) * (r2- 
2 *D* ( (D+y) /r2 ) ) + 

(gamma* 12) * ( ( D+y) /r2 ) + ( Al/r2 ) * (Pi+Pa+ ( gamma* y* 12) ) ; 

if G >= 0 % Testing for cavitation 
data_l = [data_l; G x — ( y ) ] ; 
break, end 

end 

end 

%Truncating cavitation boundaries at intersection 

index_u = f ind (data_u ( : , 3 ) <data_l ( : , 3 ) ) ; 
index_l = f ind (data_l ( : , 3 ) >data_u ( : , 3 ) ) ; 
data_u (index_u, :)=[]; 
data_l (index_l, :)=[]; 

%Plotting cavitation boundaries 

plot (data_l ( : , 2) , data_l ( : , 3) , data_u ( : , 2) , data_u ( : , 3) ) ; grid; 
axis ( [0,550,-50,5] ) ; 

title ([ 1 Cavitation Zone for a *,W_text, T lb HBX-1 Charge at a Depth of 

1 , D_text , 1 feet * ] ) ; 

xlabel ( 1 Feet 1 ) ; ylabel ( 1 Feet 1 ) ; 
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APPENDIX B. USA/LS-DYNA INPUT DECKS 



This section of this appendix provides example USA input decks for each of the 
three USA modules: FLUMAS, AUGMAT, and TIMINT. Reference 18 provides 
information concerning the various input deck variables. 

An example LS-DYNA KEYWORD input deck is included also from the full fluid 
size model. Only the first line or two of each card is included. Reference 16 provides 
information as to the meaning of each field on the cards. 



USA INPUT DECKS: 



FLUMAS DATA FOR FSP 
flunam geonam strnam daanam 


$ 


FLUNAM GEONAM GRDNAM 


DAANAM 
F F F T 


$ 


PRTGMT PRTTRN PRTAMF 


CALCAM 
F F T F 


$ 


EIGMAF TWODIM HAFMOD 


QUAMOD 
F F T F 


$ 


PCHCDS NASTAM STOMAS 


STOINV 
F F F T 


$ 


FRWTFL FRWTGE FRWTGR 


FRESUR 
F T F F 


$ 


RENUMB STOGMT ROTGEO 


ROTQUA 
F F F F 


$ 


PRTCOE STRMAS SPHERE 


ROTSYM 
F F F F 


$ 


OCTMOD CAVFLU FRWTFV 


INTCAV 
F T 


$ 


BOTREF MASREF 


0 65808 0 6186 


$ 


NSTRC NSTRF NGEN NGENF 


0 0 0 


$ 


NBRA NCYL NCAV 


0.935E-4 59428.8 


$ 


RHO CEE 


1. 1. 


$ 


CQ ( 1 ) CQ (2 ) 


o 

o 

i — i 

i 

o 


$ 


DHALF CXHF CYHF CZHF 


36. 0. 1. 0. 


$ 


DEPTH CXFS CYFS CZFS 


14.7 386.4 


$ 


PATM GRAVAC 


1 


$ 


NSRADI 


0. 0. 1 6186 1 


$ 


RADI RAD2 JBEG JEND JINC 


0 


$ 


NSORDR 


AUGMAT DATA FOR FSP 
strnam flunam geonam prenam 


$ 


STRNAM FLUNAM GEONAM 


PRENAM 
F F F F 


$ 


FRWTGE FRWTST FRWTFL 


LUMPFM 
F F F T 


$ 


FLUSKY DAAFRM SYMCON 



DOFTAB 
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F F F F 
PRTAUG 


$ 


PRTGMT PRTTRN PRTSTF 


F F F F 
CFADYN 


$ 


MODTRN STRLCL INTWAT 


11 


$ 


NTYPDA 


65808 197424 3 3 


$ 


NSTR NSFR NFRE NFTR1 


1 


$ 


NSETLC 


0 1 6186 1 


$ 


NDICOS JSTART JSTOP JINC 



TIMINT DATA FOR FSP 
prenam posnam 
resnam 
F T F F 
XXXXXX 
1 

0.0 0.00002 
T F F F 
PACKET 
F T F F 
VELINP 
F F F F 
XXXXXX 
1 

0 . 

0. -252.0 -336.0 

0. -108.0 -336.0 
201 

1 . 0 . 

0. 682E-5 
1 

60. 12.0 24.0 
99999 99999 
0 0 0 0 
NS TART 
F F F F 
XXXXXX 

0.0 -14.4 -108.0 
288. 0. 1. 0. 

14.7 386.4 
F 



$ PRENAM POSNAM 
$ RESNAM WRTNAM 
$ REFS EC FLUMEM XXXXXX 

$ NTINT 

$ STRTIM DELTIM 
$ EXPWAV SPLINE VARLIN 

$ HYPERB EXPLOS DOUBDC 

$ BUBPUL SHKBUB XXXXXX 

$ NCHARG 
$ HYDPRE 
$ XC YC ZC 
$ SX SY SZ 
$ JPHIST 
$ PNORM DETIM 
$ DTHIST 
$ CHGTYP 

$ WEIGHT SLANT CHGDEP 
$ NSAVER NRESET 
$ LOCBEG LOCRES LOCWRT 

$ FORWRT STBDA2 ASCWRT 

$ XV YV ZV 

$ DEPTH CXFS CYFS CZFS 
$ PATM GRAVAC 
$ DISPLA 



LS-DYNA INPUT DECK: 

* KEYWORD 
*TITLE 

HALF-MODEL OF FSP W/MAT 90 FLUID (ORIGINAL SIZE MODEL) 
* CON TROL_T E RM I NAT I ON 
0.080, 0, 0, 0, 0 
*CONTROL__TIMESTEP 
0 . 00002 , 0 . 8 , 0 , 0 , 0 , 1,0 
*CONTROL_OUTPUT 
1 , 0 , 1 , 0 , 0 , 1,1 
* CONTROL PARALLEL 
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«/> < J > «/> 



2 , 0 , 1 

*DEFINE_CURVE 

1 

0 , 0.00002 

0.080,0.00002 

* DAT ABAS E_H I S TORY_NODE 

225, 7 28, 157 6, 1615, 24 45, 24 71, 1005,1609 
18 68, 3815, 68 90, 19819, 22885 

* DATABASE_NODOUT 

0.00002 

* DATABASE_B I NARY_D3 PLOT 

0.0002 

* DAT ABAS E_B I N ARY_D 3 T H DT 
0.0002 

* DAT ABAS E_EX TENT_B I NARY 
0,0,3, 1,1,1, 1,1 " 

0 , 0 , 0 , 0 , 0, 0 
*BOUNDARY_USA_SURFACE 
2 , 1,0 

*INIT I AL_ DE TON AT I ON 
-1,0.0,-252.0,-336.0,0.0 

6204 . 63. 0. 28 91E-03, 0.0,-108.0,-336. 0,25072 
NODES 

*NODE 

1. 5. 99997997. 72 . 0000000. 90.0000000. 0.0 



*PART 

1 

1 , 1 , 1 , 0,1 
*SECTION_SHELL 
1 , 2 , 1 . 0, 3 

1 . 000E+00 , 1 . 000E+00 , 1.000E+00, 1.000E+00 
* ELEMENT_SHELL__TH I CKNESS 
1, 1,1005, 1121, 1122,1111 
1.000E+00, 1.000E+00, 1. 000E+00, 1.000E+00 



*PART 

2 

2 , 2 , 2 , 0,2 

*SECTION_SHELL 

2 , 2 , 1 . 0, 3 

1 . OOOE+OO , 1. 000E+00, 1 . 000E+00 , 1 . 000E+00 
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*■ -(/> -co -tr> <r> </> 



*ELEMENT_SHELL_THICKNESS 

1161, 2, 1628, 1630, 1785, 1783 

1.000E+00, 1 . 000E+00 , 1.000E+00, 1.000E+00 



* PART 
3 

3. 3. 3. 0. 3 
*SECTION_SHELL 

3. 2.1. 0. 3 

5.000E-01, 5.000E-01,5.000E-01, 5.000E-01 

*ELEMENT_SHELL_THICKNESS 

1465,3, 1576, 1712, 1734,1580 

5. 000 E- 01, 5. 000 E- 01, 5. 000 E- 01, 5 . 000E-01 



*PART 

4 

4. 4. 4. 0. 4 
*SECTION_SHELL 

4.2. 1.0. 3 

0.50, 0.50,0.50, 0.50 
*ELEMENT__SHELL_THICKNESS 
2169, 4, 2463,2466, 2300,2460 
5.000E-01, 5.000E-01, 5.000E-01, 5.000E-01 



*PART 

90 

90, 90, 90 

*SECTION_SOLID 

90,8 

*ELEMENT_SOLID 

1, 90, 310, 316, 204,305, 3670, 3671,3672, 3673 



MATERIAL CARDS 

DEFINITION OF MATERIAL 1 PSHELL bottom_si.de 

MAT_ELASTIC 

1,7. 320E-04 , 3 . 000E+07, 0.300 
* HOURGLASS 
1,5 
$ 

$ DEFINITION OF MATERIAL 2 PSHELL dryside 
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$ 

*MAT_ELASTIC 

2, 7.320E-04, 3. 000E+07, 0.300 
* HOURGLASS 

2.5 

$ 

$ DEFINITION OF MATERIAL 3 PSHELL inner_bottom 

$ 

*MAT_ELASTIC 

3,8. 200E- 04,3. 000E+ 07, 0.300 
*HOURGLASS 

3.5 
$ 

$ DEFINITION OF MATERIAL 4 PSEHLL stiffs 

$ 

* MAT_E LAS TIC 

4,8.800E-04,3. 000E+07 , 0.300 
* HOURGLASS 

4.5 

$ 

$ DEFINITION OF MATERIAL 90 PSOLID water 

$ 

*MAT_ACOUSTIC 

90,9. 35E-05, 59428. 8, 0.5, 1.0, 14. 7, 386. 4 
0.0,36.0,-96.0,0.0,1.0,0.0 
$ 

$ Face set daa_bdry 
$ 

*SET_SEGMENT 

2 

6729, 6730, 6731, 6732 



$ 

$ RAYLEIGH DAMPING 
$ 

*DAMPING_PART_STIFFNESS 

1.0. 5E-05 

2.0. 5E-05 

3.0. 5E-05 

4.0. 5E-05 
*END 
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APPENDIX C. HELPFUL FEATURES IN MSC/PATRAN 



MSC/PATRAN [Ref. 14] is a powerful finite element modeling and visualization 
tool. The program has many useful features to aid in visualizing and manipulating a finite 
element model. Different input and output formats are also supported, such as LS-DYNA 
keyword format. Basic familiarity with MSC/PATRAN is assumed. Important menu 
selections are featured in bold. 

1 . Creating GROUPS is useful for visualization of a complex mesh. The model can 
in effect be "sliced" into different sections for example and each put in a different 
group. The CREATE option is used to form a group. A group name must be 
entered. Elements to add to a group can then be selected with the mouse from the 
viewport. When selecting elements to add to a group, visible elements only can be 
selected by toggling the visible only button on the top left of the selection tool bar. 
This is only effective when using the hide view of the mesh (vice wireframe). The 
visible elements only selection feature is most useful when defining the wetted 
surface of a model. The MODIFY option provides for additions or removals from 
the target group. The target group can be changed with the CHANGE TARGET 
button. The POST option allows individual groups to be displayed in the current 
viewport. More than one group can be selected for display. The shift key must be 
held down to make multiple selections. This feature works also when selecting 
individual elements from the viewport. An individual element is selected, then, if 
additional elements are desired to be selected, they can be added to the list by 
holding down the shift key while highlighting them with the mouse. 

2. A model can be moved in set increments using the TRANSFORMATIONS 
options found under the VIEWING menu. This option provides for rotation of the 
model about one of the three axes in a set fashion. The model can also be moved in 
set increments in any of the six main directions (x,y, and z). This provides for 
precise control over the model's positioning and aspect. 

3. Multiple viewports can be created and posted using the VIEWPORT menu. The 
CREATE option allows the user to create and name a new viewport. Each created 
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viewport can be posted/unposted using the POST menu. The MODIFY menu 
provides for the current and default viewports to be changed. The current viewport 
is the "active" viewing window and it is where all actions performed will take effect. 
Groups can also be posted/unposted from this menu option. This effects the current 
viewport. The TILE option automatically places two displayed viewports side by 
side. 

4. The LIST creation option under the TOOLS menu provides the means to find 
objects (elements, nodes, etc.) with a common ATTRIBUTE or ASSOCIATION 
with other objects (such as groups or certain elements). FEM or geometry can be 
chosen. Nodes or different element types can be specified for association or attribute 
selection. For example, an empty group can be created (the group name is entered in 
the group create dialog box and apply is depressed; no elements are selected from 
the viewport) and then a list can be created of all the elements with a common 
attribute, such as a particular material property or property set. The list of elements 
can then be added to the newly created group and subsequently displayed. The list 
can also be added to any existing groups. The list can be output to one of two sub- 
windows, A or B. The destination window is selected via the radial buttons at the 
bottom of the create list window. The list tool can be a very powerful asset. 

5. The DISPLAY menu option provides a number of menu selections. The 
PLOT/ERASE option provides the means to "unclutter" the display viewport. 
Specific element types can be selected and then erased from the display. The erased 
elements are not deleted from the model; they are only removed from view. They 
can then be re-posted to the viewport by selecting the PLOT ALL POSTED FEM 
OPTION. All the objects in the viewport can also have their labeling toggled on 
and off via the appropriate menu selection. The object color can also be changed 
from the set default color. This is accomplished by simply clicking the small color 
patch next to the object name in the appropriate label selection menu. A small 
window of color choices will pop-up and a new color can be selected and applied 
the object. 
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6. The FINITE ELEMENT radial button has a number of useful mesh creation and 
diagnostic tools, one of which is the VERIFY option. By selecting ELEMENTS 
and NORMALS, element normal vectors can be displayed and even reversed if 
required. For the reverse option, a reference element must be designated. The 
ELEMENTS/DUPLICATES option allows any duplicate elements to be 
highlighted and deleted if the user desires. Either the higher or lower ID number 
element can be selected for deletion. EQUIVALENCE allows a tolerance to be set 
and any nodes falling within the tolerance will be merged and the database 
numbering reset. Nodes can be designated for exclusion from the equivalencing. 
The RENUMBER option allows nodes and/or elements to be renumber starting 
with a user specified number. The SHOW option displays the selected node's 
coordinates and ID number, and for elements the ID number, type, and property set 
is displayed. The selected nodes and elements are highlighted in the viewport if the 
ID number is input by hand rather than selected with the mouse. 

7. Selecting the LOADS/BCS radial button allows creation of pressure load to 
define the wetted surface. When inputting the pressure set data, one has choices of 
top or bottom of the element (for two-dimensional elements). The correct choice is 
top, since this is the side with the normal vector pointing outward. For three- 
dimensional elements, such a fluid elements, the pressure load can be applied to free 
faces only (a button on the selection tool bar). This is useful for defining the outside 
surface of the fluid mesh for the DAA boundary. The visible only button should also 
be depressed when defining a wetted surface. 

8. Results from a NASTRAN analysis can be input into PATRAN by selecting the 
ANALYSIS radial button, followed by choosing READ OUTPUT2 under the 
action menu. The desired results file name can then be selected and read in by 
clicking the apply button. 

9. The imported results can then be viewed using the RESULTS button. If the 
results are in a form that can be plotted in an xy-plot, then the results type should be 
changed from BASIC to ADVANCED. Under RESULT CASE OPTIONS, the 
desired results should be selected (highlighted), then the GET RESULTS bar 
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depressed. Under the PLOT TYPE menu, XY-PLOT can be selected. PLOT 
TYPE OPTIONS is used next to assign the global variable (usually time) and then 
the desired y-variable (such as displacement, velocity, and acceleration, depending 
on the results read in) can be selected to be plotted. The nodes whose response is 
desired can be input. Once the results curves are plotted, they can be further 
manipulated using the XY data radial button. 

10. The XY data button has the normal menu selections, such as CREATE, 
DELETE, MODIFY, POST. These choices have options for both xy-windows and 
curves. Using create, new xy-windows can be generated and displayed with the post 
option. Curves cannot be generated here, they must be generated as in step 9 above. 
Curves can be deleted here and posted/unposted to/from different xy-windows. 

1 1 . ASCII files can be created from the generated curve data. Under the xy menu, 
MODIFY/CURVE should be selected. The desired curve is then selected 
(highlighted) from the middle window (which displays all of the curves currently 
created in the database file). The DATA FROM KEYBOARD option should be 
chosen and the WRITE XY DATA TO FILE button should be clicked on. Once 
apply is clicked, a filename can be input for the ASCII file and the destination 
directory selected. PATRAN does put a small text header in the file above the first 
xy-data pair. This header must be deleted prior to importing the file into UERD 
Tool. 

12. The PRINT option can be found under the FILE menu. The user must be 
careful when printing to select what is to be printed. The top bar of the print menu 
can be toggled between either CURRENT VIEWPORT printing or CURRENT 
XY- WINDOW printing. ALL VIEWPORTS or ALL XY-WINDOWS can also 
be chosen. Color or black&white can be toggled under the print OPTIONS menu. 
Under options is where the print to file option is chosen and a filename is input. 
PATRAN can output postscript or encapsulated postscript files. The PAGE SETUP 
menu provides for selection of page orientation and size, as well as turning borders 
on and off. The output is not created until the apply button is clicked. It should also 
be noted that the output file name need not be changed for subsequent print views. 
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PATRAN automatically appends a sequentially numbered suffix to the postscript 
filename with each print output. 

13. NEUTRAL files can be generated by selecting the EXPORT option under the 
FILE menu. Neutral files can be imported into PATRAN using the IMPORT 
option. PATRAN by default looks for a .out extension on neutral files. 
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APPENDIX D. FLUID MODELING USING TRUEGRID 



This appendix covers the procedure for creating a fluid finite element mesh using 
TrueGrid's extrusion feature: the ELUDE command. The basics of using TrueGrid will not 
be covered here and some familiarity with the code is assumed. Additional information can 
be found in the TrueGrid user manual [Ref. 12]. 

Essentially the BLUDE feature pulls or "extrudes" the structural mesh through a 
"guide" mesh mated to the structural wetted surface in the form of a block part. The block 
part is actually attached to a surface definition created from a faceset of the wetted elements 
of the structural mesh. The resulting extruded mesh matches exactly to the structural mesh, 
a prerequisite for successful fluid modeling. 

The extrusion procedure is as follows, with important commands and menu 
selections denoted in bold and all capital letters for emphasis: 

1. A structural model must be created. TrueGrid can be used or the READMESH 
command can be used to input a mesh from another code format, such as LS-DYNA 
or NASTRAN. It is very important to remember though, that when TrueGrid reads 
in a finite element mesh from an outside code format, it renumbers ever element and 
grid point. Therefore, once the mesh is through being manipulated in TrueGrid, and 
it is written an output file, the grid point and element ID numbers will not match 
between the original and newly output model from TrueGrid (even if the original 
model was not modified in TrueGrid). 

2. The elements of the structural model that will be in contact with the fluid, i.e. the 
wetted surface, must be grouped into FACESETS. This option can be accessed 
from the environment window under the PICK option by choosing the SETS 
button. The FACES button should be selected. Each "face" of the structural model 
should be put in a separate faceset, meaning each side, bottom, bow, and stem 
should be grouped individually. The reason for this will be clear once the procedure 
is understood and used. The HIDE drawing mode vice WIREFRAME should be 
used for the mesh to ensure that only the visible elements are picked. This will make 
faceset selection must easier, since it must be done by hand using the lasso tool 
guided by the mouse. The four-node selection option is the best to use when 
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choosing the faceset. This means that four nodes of an element must be within the 
selection lasso for the element to be added to the faceset. The selected elements will 
be highlight in white. If some elements are selected that are not desired in the 
particular set, they can be easily selected and removed; using the one node selection 
option is best for this operation. The REMOVE button should be pushed also. The 
set must be named and saved once selected. 

3. The SURFACE menu SD (surface definition) option should be chosen next. A 
surface number must be input. The faceset option should be selected from the end of 
the surface options list and the name of the desired faceset should then be input. 
This step converts the named faceset into a surface definition. The new surface will 
be displayed in red in the physical window. 

4. Next, the PARTS menu should be selected and the BLUDE option chosen. 
Using this option, the user creates a block part that will be attached to the above 
created surface. This block will serve as the "guide" for the extrusion of the 
structural mesh; therefore, the block's mesh must match the structural mesh or be of 
finer quality in order to get a quality extrusion; an exact match is not required 
however. This block part is created in the same way as a block using the BLOCK 
command. The blude command requires two additional inputs, however. First, the 
face of the block where the extrusion begins must be input. This is simply the face 
closest to the structure. Next, the name of the faceset to be extruded must input. 

5. The block part created must now be attached to the surface created in step 3. It 
can be attached using any of TrueGrid's available options. The easiest being 
selection of the face to be attached and then selecting the surface and clicking the 
PROJECT button in the environment window. This will work for simple cases, but 
a complex surface may require use of other TrueGrid methods. 

6. The interface of the extrusion mesh and the structural mesh should be carefully 
examined. Orthogonality of the fluid and structural mesh is a must (next to the 
wetted surface) and should be verified; TrueGrid's DIAGNOSTICS menu provides 
the necessary tools. The block mesh can be modified as needed using various 
TrueGrid tools to ensure a quality mesh is constructed for the extrusion; two 
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examples of useful tools are the mesh relaxation algorithms and use of a cubic spline 
to added curvature to the block mesh edges. Material properties can be assigned to 
the mesh also, just as with any other part in TrueGrid. In short, the extrusion mesh 
should be treated as any other part one would create in TrueGrid; all of the same 
options are available. 

7. Once the user is satisfied with the extrusion mesh, the MERGE command 
should be used to end the PARTS phase and actually perform the extrusion. The 
result will be a fluid mesh, which matches exactly to the structural mesh. The mesh 
will consist of 8-noded solid elements. The STP option can be used also if required 
to ensure that the fluid mesh is merged with the structural mesh and there are no 
duplicate nodes. Additionally, prior to merging, the extrusion mesh can be 
replicated using the LCT and LREP commands. This will only be effective if the 
model is symmetric. Using these part replication features, the user only has to build 
one-half of the extrusion mesh. 

8. Additional extrusions can be performed, including on any newly extruded mesh 
surfaces. This must usually be done to fully form a fluid mesh around the structural 
model. 

9. Postscript images of the model and the mesh can be made using the 
POSTSCRIPT command. The command postscript is given at the command 
prompt with a the desired output filename. The DRAW button in the environment 
window should then be clicked to redraw the image. This creates the postscript file. 
Additional files will be generated as long as the command is active and the model is 
manipulated in any way so that it must be regenerated in the display window. The 
postscript command can be turned off by typing POSTSCRIPT OFF. One 
additional command that is quite useful in generating quality image files is the 
RESO command. The reso command is entered prior to the postscript command. 
The syntax is the command followed by a number that is the desired resolution. This 
is system limited. A resolution of 2300 has been used with success. The postscript 
files generated are black and white only, but they consist of vector data. This means 
the images generated are crisp and very accurate. 
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APPENDIX E. USEFUL FEATURES IN LS-TAURUS 



TAURUS is an interactive post-processor and three-dimensional visualization tool 
for LS-DYNA [Ref. 16]. TAURUS has many useful features for viewing and manipulating 
model output responses obtained from LS-DYNA/USA simulations; a few of these features 
are detailed here. All TAURUS commands are fully documented in Appendix K of 
Reference 16. Basic familiarity with TAURUS is assumed. 

1. A given state can be displayed in TAURUS by the following syntax: s (state 
□umber) frin (fringe number). Note that the parentheses are not included (this 
syntax will be used to illustrate all commands); only the appropriate number is put 
in the command. In underwater shock simulations, pressure fringes are of concern; 
fringe 8 is pressure The fringe numbers for other variables are listed in a set of 
tables in Reference 16. 

2. Animation can be easily set-up with the command: r (starting state number) 
(final state number) (step) frin (appropriate fringe code). The animation can be 
set-up to run between any desired starting and ending state. The gif command 
allows a gif image to output of the currently displayed state in the visualization 
window. TAURUS names the file by default, pict#.gif (where the # is a consecutive 
integer starting at one for each gif image created during the current TAURUS 
session). The command noborder removes the border and text information from the 
visualization window. Logo toggles the drawing of the TAURUS logo. The 
command cb can be used to change the background color from the default black. 
After entering the command, the user must input the amount of red, green, and blue 
to be mixed for the desired background color. As an example, 0,0,0 is black, and 
1,1,1 is white. Other colors can be made by experimentation. The command cline is 
used to change the mesh line color on the screen. The default is white. As with cb, 
the amount of red, green, and blue must be input. 

3. A video tape of the resulting animation can be made by selecting the video out 
application from the SGI desktop list of media tools. Note that this feature must be 
installed and a VCR and TV must be hooked up to the SGI machine. The NPS SGI 
OCTANES have this option installed. The outline of a box will be visible on the 
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screen. This is the capture window outline for what is sent to the VCR/TV. The size 
of this window can be changed between two sizes: full screen and a somewhat less 
than full screen size. The size is changed by clicking on the video out gray box that 
appears on the screen. The flicker filter can also be toggled on and off in the same 
way, and the function can be turned off here too. The full screen size does not give 
the best output; it is better to use the smaller size for better quality on the TV end. 
The VCR can then be used to record the images transmitted from the computer. 

4. Phase II of TAURUS is used to plot time histories of desired variables for 
elements and nodes. The procedure for making the plot is as followings (example is 
for element time history): first the command elem (number of elements) (element 
numbers) is used to designate the desired elements; next gather is used to read the 
element data into memory; once the gathering is complete the desired time history 
can be plotted with etime (fringe number) (number of elements to plot) (element 
numbers to be plotted). The commands aset (min) (max) and oset (min) (max) 
can be used to set the x and y axes appropriately for the desired range. The keep 
command can be used to generate an ascii file with the plotted data in it. This 
command must be issued prior to the etime command and once entered the user is 
prompted to enter a filename for the output file. The ASCII file generated is in two 
column (x-y) format and contains the entire range of data for the element fringe 
requested (even if the data is only partially displayed in the plot window by use of 
the axis setting commands). It should be noted that this ascii file does contain some 
text header information above the data columns and the string "endplot" is placed at 
the end of the data. 
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